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We derive the effective channei for a iogicai qubit protected by an arbitrary quantum error- 
correcting code, and derive the map between channeis induced by concatenation. For certain codes 
in the presence of single-bit Pauli errors, we calculate the exact threshold error probability for perfect 
fidelity in the infinite concatenation limit. We then use the control theory technique of balanced 
truncation to find low-order non- asymptotic approximations for the effective channel dynamics. 
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Quantum Error Correction Q and Fault-Tolerant 
Computation Q have demonstrated that, in principle, 
quantum computing is possible despite noise in the com- 
puting device. Analyses in these areas, as in many 
physics disciplines, rely on asymptotic limits and expan- 
sions in small parameters. However, realistic devices will 
be of large but finite scale and real parameter values 
may not be sufficiently small; thus methods valid out- 
side these limits are desirable. Here we demonstrate the 
use of model reduction ^J, a control theory technique, 
for studying large but finite quantum systems. 

We will first derive the effective channel describing an 
encoded qubit 's evolution under arbitrary error dynam- 
ics, and propagate these results to concatenation schemes 
jjj when the dynamics do not couple code blocks. The 
calculations for both finite and asymptotic concatenation 
are simple, but the resulting expressions for the effective 
channel dynamics in the finite case are cumbersome and 
ill-conditioned. We will then use model reduction to find 
low-order approximations for these dynamics. Among 
other results, we derive a model for a 9 4 -qubit system re- 
quiring only 23 degrees of freedom for accuracy ~ 10~ 3 . 

We consider the following formulation of quantum er- 
ror correction: a single-qubit state po is perfectly encoded 
in a multi-qubit register with state p, which evolves via 
some error dynamics described by a linear map p(0) <— > 
p(t) = £[p(0)). (For master equation evolution p = C[p], 



£ 



) At time t, a syndrome measurement is made 



and the appropriate recovery operator is performed; we 
assume the measurement and recovery are noiseless. This 
process returns the system to the codespace |(| , and thus 
its state can be described by a single-qubit state pf. 
Denote the probability- weighted average of pj over syn- 
drome measurement outcomes by ~pj. For a given code 
and error model, we wish to know ~pj(t): we may then 
compare po and ~pj(t) via a desired fidelity measure. 

The logical qubit po may be parameterized by the 
expectation values (I)o, (X) Q , (Y) and (Z) , with 
{/, X, Y, Z} the usual Pauli matrices. (Of course (J) = 1, 
but it will be convenient to include this term.) The log- 
ical qubit po is encoded by preparing the register in the 
state p(0) = (J) £j + (X) E X + (Y) E Y + (Z) E Z where 
E a acts as \u on the two-dimensional codespace and van- 
ishes elsewhere. These operators are easily constructed 



from the codewords: for example, given the encoding 
|0> ^ |0), |l) ^ |T), we have E x = |(|0)(T| + |T)(0|). 
(For a stabilizer code [[[) with stabilizer S — {Si} and 
logical operators a, the codespace projector is Pc — 
■4^5, and E a = \P c cr.) E.g., for the bitffip code 
m given by |0) i-> |000), |1) ^ |111), 



Ei = |( III + IZZ + ZIZ + ZZI ) 
E x = I( XXX - XYY - YXY - YYX ) 
E Y = I( -YYY + YXX + XYX + XXY ) 
E z = |( ZZZ + ZII + IZI + IIZ ). 



(1) 



For the trivial code (po "encoded" as itself in a single- 
qubit register) we have simply E„ = i<r. 

After the action of £, recovery yields the expected log- 
ical state ~p]. As with po, we parameterize ~p] by the 
expectation values (I) j, (X) j, (Y) j and (Z) j. The (a) j 
may be written as expectation values of operators on the 
register prior to recovery: letting {Pj} be the syndrome 
measurement projectors and {Rj} be the recovery oper- 
ators, {a) j = ti(D a p) where D a = 2 V, I'" if E 7 li , I',. 
(For a stabilizer code, Rj is some Pauli operator of low- 
est weight leading to syndrome measurement Pj, and 
D <y = J2i fiaSiV, with f ia = J2jV(Si,Rj)v(Rj,v) 
where q) = ±1 for pq = iqp.) For the bitffip code, 
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For the trivial code we have simply D a = a. 

To describe the evolution of the encoded logical bit, 
we compute the effective channel Q taking po to ~pj. 
Q may be written as the linear mapping of po = 
«/) , (X) , (Y) , (Z) ) to pj = ((!}?, (X) f , (Y) f , (Z) f ): 



(ff)/ = tr D a £ 



(3) 



If £ is completely positive |Q, it follows that Q is as well. 
To obtain a matrix representation of Q, let 



tv{D a £[E a ,]). 



(4) 
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Then pj = Gpo, and the fidelity of a pure logical qubit 
under this process is ^p^Gpa- Note that the dynamics £ 
need not be those against which the code protects. 

Now consider concatenated codes jj]]. In the concate- 
nation of two codes, a qubit is encoded using the outer 
code C ont and then each of the resulting qubits is en- 
coded using the inner code C m . Though not necessar- 
ily optimal, a simple error-correction scheme coherently 
corrects each of the inner code blocks, and then corrects 
the entire register based on the outer code. We denote 
the concatenated code (with this correction scheme) by 
C out (C in ). Given the effective channel Q describing C in 
under some dynamics, and a desired C out , we now con- 
struct Q describing the evolution under C out (C m ). 

Let C out be an iV-bit code and C in be an M-bit code. 
We assume that each M-bit block evolves according to 
the original dynamics £ and no cross-block correlations 
are introduced, thus the evolution operator is £ = £ ® 
£ ® . . . ® £. Each M-bit block represents a single logical 
qubit encoded in C m ; as the block has dynamics £, this 
logical qubit's evolution is described by Q. Therefore 
the effective evolution operator for the logical bits in the 
codeword of C out is £ cB = Q ® Q ® . . . ® Q. 

Operators on N qubits may be written as sums of ten- 
sor products of N Pauli matrices; we may therefore write 
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D': 



{I,X,Y,Z} 



(M ( 5 ) 



(6) 



(For stabilizer codes, the a and (3 coefhcients are easily 
found.) Substituting £ eS , and L>° ut into (|) and 

noting that tr(i/j-£?[g#j]) = Gv jH yields 



Goo' 





(7) 



Thus the matrix elements of G can be expressed as poly- 
nomials of the matrix elements of G, with the polynomial 
coefficients /^.iQ^.-j. depending only on the E a and D a 
of the outer code. For a given outer code C, denote the 
concatenation map G > G by il c . It can be shown that 
f2 preserves complete positivity. 

More generally, the inner process G can represent any 
linear evolution of a logical qubit, and thus we may speak 
of concatenating a qubit process with a code. E.g., if 
G describes some qubit dynamics, then f2 (£7) describes 
the effective dynamics of encoding by C with the un- 
corrected dynamics G acting on each register bit. This 
method only requires that the outer code's logical qubits 
be decoupled. (Above we assumed £ e ff = G <8> ■ ■ ■ <8> G', for 

£ c ft = G {1) <8> . . . <8> G {N \ replace Q Vjilj with Gi% 3 in (§.) 

We may characterize both the finite and asymptotic 
behavior of any concatenation scheme involving the codes 



{Cfc} by computing the maps Qr k . Then the finite con- 
catenation scheme C±(C2(- ■ ■ C n . . . )) is characterized by 
Qd(c a (...c n ... )) = Qd ( d(, . .n c » . . . )). We expect the 

typical fl c to be sufficiently well-behaved that standard 
dynamical systems methods || will yield the t — > oo limit 
of (Q C Y; one need not compose the (fl c 'Y explicitly. 

We now consider certain concatenation schemes when 
the symmetric depolarizing channel jlj acts on each reg- 
ister qubit. This channel is described by G dcp (t) diagonal 
with entries (1, e -7 *, e _7t , e~ 7 *). From trace preservation 
Gn is always 1, so let [x, y, z] denote G diagonal with en- 



[e- 7 *,e" 7t , 



-7*1 



tries (l,x,y,z); then G dep (t) 

Suppose more generally we are given a qubit process 
described by G — [x,y,z], and wish to concatenate this 
process with the bitflip (bf) code. Using fflt) and the cod- 
ing operators (|l|) and (0), we find that fl M (G) describing 
the concatenated evolution is also diagonal: 



Q bt ([x,y,z}) 



[x 3 , \x 2 y 



(8) 



((H) could also be found by using the Heisenberg picture 
to evaluate (|).) Writing |±) = -^(|0)±|1)), the map fiP f 
for the phaseflip code |±) t—> |±±±) Q is similar. (Any 
stabilizer code Cs preserves diagonality: f2 s ([x, y, z]) is 



Goo' — 3oo' <gi ^ t fioX 



wx(SiX) y wy(SiY) z wz(SiZ) ^ 



with w a (p) the a- weight of a Pauli operator p (e.g. 
wx(XYX) = 2) and fi a as previously defined.) 

The concatenation phaseffip(bitffip) yields the encod- 
ing |±> i — ► -ig(|000)±|lll)) <33 , which is the Shor nine-bit 

code @. Thus Q sho1 = fiv f (tt hf ), and Q sho1 ([x , y , z}) = 
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-Hf^-i, 3 ) 3 , (lz-^) 



(10) 



Now consider the Shor code concatenated with itself 
£ times, and let |s<(f),&(t),2<(t)] = (ft Shol f {G dcp (t)). 
The functions ae(t) approach step functions in the limit 
I — » oo (e.g., see Fig. Q); denote these step functions' 
times of discontinuity by i*. Thus in the infinite con- 
catenation limit, the code will perfectly protect the (a) 
component of the logical qubit if correction is performed 
prior to t*. We call t* the a- storage threshold. 

We calculate the t* by finding the I — » oo limit oiae.it). 
Writing (JToj) as [Qi(x), Qz(x,y, z),Q 3 (z)}, the map z ^ 
Q^{z) has stable fixed points at and 1 and one unstable 
fixed point z* on (0, 1); numerically solving Qz{z) = z 
yields z* w 0.730. The plots of Zg(t) all intersect at 
z(t^) = z*, and the step function limit follows from the 
stability of and 1. Inverting e -7 * = z* yields t* z . Q\ has 
similar features with x* w 0.900, and a similar analysis 
yields t* x . For this code, one can show ty — niin(tjt, tjj). 

We may also phrase the thresholds in the language 
of finitely probable errors. The expected evolution of a 
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FIG. 1: zt (t) for Shor £ concatenation under the depolarizing 
channel; the fidelity of an encoded Z eigenstate is i(l + 2jp(i)). 
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0.1050 


0.3151 


0.1618 


0.2150 


0.1383 


0.2027 


pth 


0.0748 


0.1121 


0.0969 


0.1376 



TABLE I: Code storage thresholds, 
qubit subjected to a random Pauli error with probability 



pisp^ (l~p)p+^XpX + ^YpY+^ZpZ. This channel 
is described by g Pauli (p) ='[1 - fp, 1 - §p, 1 - |j>]. As 
gPauiid^ _ e - T t^ _ gdep^^ in the innnite concatena- 
tion limit with G Flluh (p) acting on each register qubit, the 
logical qubit 's (tr) component will be perfectly protected 
if p < p* = 1(1 — e f*"-). Define the threshold prob- 
ability pth = min{p*}; for p < p t h, all encoded qubits 
are perfectly protected in the infinite concatenation limit. 
Values for t* and p t h appear in Table |. 

For comparison, we derived thresholds for three other 
codes. Another version of the Shor code is given by |0) i— > 
^(|000) + |111))® 3 , |1) i ^ ^g(|000> - |111))® 3 ; call this 

code^Shor'. Let [x' t (t), y[{t), 2j(t)] = (Q sW )* (0 de P(i)). 
The y' e (t) approach a step function as t — > oo, but x^(t) 
and approach a limit cycle of period 2, interchanging 
step functions with different discontinuities at every iter- 
ation of n ShoT . Considering instead the limit of iterating 
(O shor ) 2 permits an analysis as for fi shor . The Steane 
seven- bit code may be treated similarly to the Shor 
code, and the symmetries of the Five-Bit code Q lead to 
a simple analysis. Results are summarized in Table [jj 

We now return to the Shor code under the depolariz- 
ing channel, and consider the finite concatenations de- 
scribed by the functions o^(t). The ai(t) have the form 
X)i bie~ ant with the a, positive integers and the bi ratio- 
nals. For 7 = no errors occur, thus J2i bi = 1. 

Explicit calculation of the at has several disadvantages. 
First, the number of terms in these series grows approx- 
imately as 9 e (see Table ||(a)). Though not nearly as 
severe as for the number of elements in the full-system 
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TABLE II: (a) Terms in exact series for ai(t). (b) Order of 
iteratively reduced realizations for ae(t). 

density matrix (2 2 ' 9 ), this growth is still too rapid to be 
practical. Only a small portion of the terms in these se- 
ries have \bi\ < 1, thus one cannot meaningfully truncate 
the series without introducing significant error. 

More seriously, the magnitude of the bi grows rapidly: 
e.g., \bi\ > 10 60 for 65 of the 352 terms in z 3 , and double- 
floating point precision no longer yields Y]- bj = 1. To 
efficiently generate plots of the 3^(t) we repeatedly apply 



Shor 



to numerical values of 



-7* 



-7* 



-7*1 



for all de- 



sired times t. However, this leaves us without a dynamic 
model for the evolution of ~p](t) ■ 

Given Q(t) — [x(t),y(t),z(t)], for each a(t) we will 
seek a square matrix A ai column vector B a and row 
vector C a such that a(t) ~ C a e Aat B a . For n x n A a , 
we say (A a , B a , C' a ) is an order n realization of a(t). 
(These methods may be generalized to non-diagonal Q{t) 
by seeking matrices A, B and C of sizes n x n, n x 4 
and 4 x n respectively such that Q(t) w Ce At B.) For 
a(t) = J2i bie~ ant , we can exactly realize a(t) by choos- 
ing A a diagonal with entries —0^7, B a with entries bi, 
and C a = (1,1,..., 1). If the ai are distinct and the bi 
non-zero, this realization is minimal: there is no lower- 
order exact realization of a(t). 

To find approximate lower-order realizations we use the 
model reduction technique of balanced truncation H Q. 
Consider a system with time-varying input u(t) G M, 
state x(t) G M™, dynamics x = Ax + Bu, and output 
y(t) = Cx G R; if u = S(t), y = Ce At B for t > 0. Note 
that (A,B,C) -> (TAT~ 1 , TB , CT - 1 ) leaves the map 
VP : u(t) 1 — * y(i) unchanged. An arbitrary truncation of 



an aiz 
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Cl C 2 



state-space dimensions, e.g. 

([an], [61], [ci]), may yield a radically different map \E'. 
However, we may numerically construct a balancing 
transformation T such that in the balanced system, a 
non-negative real Hankel Singular Value (HSV) hi is as- 
sociated with each dimension of the state-space M™ . Re- 
moving all dimensions with hi = yields a minimal real- 
ization; further truncating dimensions with small HSVs 
introduces a small error in 'J which, in an appropriate 
norm, is bounded by the sum of the truncated HSVs @. 

Writing the series for a(t) as minimal realizations, we 
can balance and calculate their HSVs. In Fig. || we see 
the HSVs for z(t) after each level of bitflip and phaseflip 
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FIG. 2: Largest HSVs for exact realization of z(t) at levels of 
3-qubit concatenation (17 smaller values for Shor 2 not shown). 
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FIG. 3: Exact 22(f), and approximations that result from 
balanced truncation. The order 4 approximation is only dis- 
tinguishable from the exact function on the inset. 

concatenation up to pf(bf(pf(bf))) = Shor 2 . Note that 
the number of non-zero HSVs grows rapidly at each level 
of concatenation, but the number of HSVs above any 
hmin grows slowly. (x(t) and y(t) give similar results.) 

Consider 2a (t), with minimal realization of order 37: 
the first five HSVs are (2.5 x lO" 1 )/^, (3.7 x Kr 2 )/ 7 , 
(5.3 x 10~ 3 )/7, (6.0 x lCT 4 )/7, and (5.4 x 10^ 5 )/ 7 . Trun- 
cating all but the four most significant dimensions yields 
an approximation almost indistinguishable from the ex- 
act 22(f); truncating further to realizations of order 3 and 
2 only mildly degrades the approximation (see Fig. |J). 




FIG. 4: Approximation error for ze(t) generated by iterative 
reduction with h mln = (4 x l(T 5 )/7. 



Given realizations for the a(t), we may construct real- 
izations for polynomials of the cr(t) as follows. Given 
/(f) = Cf-e^Bf and g(t) = 
tion f(t)g(t) is realized by (Af 



C g e 



B n , the func- 



C g ). The function /(f) 



B f 
B„ 



t + 1 ® Ag,Bf ® 

- 17(f) is realized by 

[Cf C g ]j. For a scalar a, the function 

a f{t) is realized by [Af, Bf, aCf). Composing these op- 
erations allows any polynomial of the a(t) to be realized, 
and thus we may directly apply the f2 to realizations. 

For £ > 2 it is impractical to construct the exact 
ae(t) and then apply balanced truncation. Instead, we 
build approximate realizations for the ag(t) using an it- 
erative approach. Begin with minimal realizations for the 
a (t) describing £ dcp (f). Alternately apply Q bf and f2 pf 
to these realizations; after each concatenation, balance 
and truncate dimensions with HSVs less than some /i m i n . 
Choosing /i m i n = (4 x 10~~ 5 )/j yields realizations with 
orders shown in Table ||(b). Comparing to Table [n](a), 
we see the resulting order reduction is dramatic. 

Fig. ^ shows the differences between the exact Ze(t) 
and the results of the iterative reduction method. Re- 
sults for approximating xt(t) and yt{t) are similar. Up 
to eight 3-qubit concatenations, the worst errors |AS^(t)| 
are only « 3 x 10 -3 . Note that the errors appear to have 
characteristic frequencies; the error is analogous to the 
ringing in frequency- limited approximations of step func- 
tions. To good accuracy the mutual intersection points of 
the xe(t) and of the zt{t) are preserved; this is expected 
as the concatenation polynomials are unchanged. 

These results suggest balanced truncation is a powerful 
approximation tool in quantum settings. Future work 
will further investigate the iterative reduction method, 
and attempt to find bounds on the approximation errors. 
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We may choose such a recovery procedure even if the code 
is imperfect. 

Balanced truncation for u = 5{i) will be discussed else- 
where. 



